library(tidyverse)
library(sensemakr)
library(estimatr)
library(texreg)
library(corrr)
library(lmerTest)
library(lme4)
library(interplot)
library(margins)
library(stargazer)
library(patchwork)
library(countrycode)
library(ggrepel)
library(regweight)

## Functions

zero1 <- function(x, minx=NA, maxx=NA){
  res <- NA
  if(is.na(minx)) res <- (x - min(x,na.rm=T))/(max(x,na.rm=T) -min(x,na.rm=T))
  if(!is.na(minx)) res <- (x - minx)/(maxx -minx)
  res
}



load("cleaned_data/Study3_data.RData")

ggplot(countrydata,aes(x=kinship_score,y=zero1(scores)))+geom_smooth()+theme_bw()+ylab("LGBT Liberalism")+xlab("Kinship Tightness Index")+geom_text(aes(label=isocode),position = position_jitter(height=.1,width=.1))+ylim(0,1)+xlim(0,1)
ggsave("output/lgbt.png",width=8,height = 1, dpi = 600, device = 'png')

a1 <- (lm(zero1(scores)~kinship_score+continents_africa+continents_asia+continents_europe+continents_oceania+continents_north_america,countrydata))

a2 <- (lm(zero1(scores)~kinship_score+continents_africa+continents_asia+continents_europe+continents_oceania+continents_north_america+ln_time_obs_ea+small_scale+agriculture+malariaindex,countrydata))

a3 <- (lm(zero1(scores)~kinship_score+continents_africa+continents_asia+continents_europe+continents_oceania+continents_north_america+ln_time_obs_ea+small_scale+agriculture+malariaindex+hdi+cath1900+muslim1900+malariaindex + urbfrac1900,countrydata))


b1 <- (lm((total)~kinship_score+continents_africa+continents_asia+continents_europe+continents_oceania+continents_north_america,countrydata))

b2 <- (lm((total)~kinship_score+continents_africa+continents_asia+continents_europe+continents_oceania+continents_north_america+ln_time_obs_ea+small_scale+agriculture+malariaindex,countrydata))

b3 <- (lm((total)~kinship_score+continents_africa+continents_asia+continents_europe+continents_oceania+continents_north_america+ln_time_obs_ea+small_scale+agriculture+malariaindex+hdi+cath1900+muslim1900+malariaindex + urbfrac1900,countrydata))

stargazer(a1, a2, a3, b1, b2, b3 ,keep = "(kinship_score)",add.lines = list(c("Continent Fixed Effects","Yes","Yes","Yes","Yes", "Yes", "Yes"),c("Control for Ethnic Group","Yes","Yes","Yes","Yes","Yes","Yes"),c("GDP/HDI/Religion Controls","No","No", "Yes","No","No", "Yes")),dep.var.labels = c("Anti-LGBT Laws","Social Expenditures"),keep.stat = c("rsq","n" ),out = "output/lgbtlaws1.tex",no.space = T,covariate.labels = c("Ancestral Kinship Strength"))

ggplot(countrydata,aes(x=kinship_score)) + geom_histogram() + theme_bw() + xlab("Kinship Tightness")
ggsave("output/histogram_lgbt.png", dpi = 600, device = 'png')

a2 <- (lm(zero1(scores)~kinship_score+continents_africa+continents_asia+continents_europe+continents_oceania+continents_north_america+ln_time_obs_ea+small_scale+agriculture+malariaindex+muslim1900,countrydata))

sensitivity1 <- sensemakr(a2, treatment = "kinship_score",
                         benchmark_covariates = "muslim1900",
                         kd = 1:6)
summary(sensitivity1)
plot(sensitivity1)
robustness_value(a1, covariates = "kinship_score")

a22 <- (lm(zero1(scores)~kinship_score+continents_africa+continents_asia+continents_europe+continents_oceania+continents_north_america+ln_time_obs_ea+small_scale+agriculture+malariaindex+muslim1900+hdi,countrydata))

sensitivity2 <- sensemakr(a22, treatment = "kinship_score",
                          benchmark_covariates = "hdi",
                          kd = 1:10)
summary(sensitivity2)
plot(sensitivity2)

sensitivity3 <- sensemakr(b3, treatment = "kinship_score",
                          benchmark_covariates = "muslim1900",
                          kd = 1:6)
summary(sensitivity3)
plot(sensitivity3)

sensitivity4 <- sensemakr(b3, treatment = "kinship_score",
                          benchmark_covariates = "hdi",
                          kd = 1:35)
summary(sensitivity4)
plot(sensitivity4)

